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Abstract 

Metamaterials are artificial composite structures designed for controlling waves or fields, and exhibit interaction phenomena 
that are unexpected on the basis of their chemical constituents. These phenomena are encoded in effective material parameters 
that can be electronic, magnetic, acoustic, or elastic, and must adequately represent the wave interaction behaviour in the 
composite within desired frequency ranges. In some cases for example, the low frequency regime - there exist various 
efficient ways by which effective material parameters for wave propagation in metamaterials may be found. However, the 
general problem of predicting frequency-dependent dynamic effective constants has remained unsolved. Here, we obtain novel 
mathematical expressions for the effective parameters of two-dimensional metamaterial systems valid at higher frequencies and 
wavelengths than previously possible. By way of an example, random configurations of cylindrical scatterers are considered, in 
various physical contexts: sound waves in a compressible fluid, anti-plane elastic waves, and electromagnetic waves. Our results 
point towards a paradigm shift in our understanding of these effective properties, and metamaterial designs with functionalities 
beyond the low-frequency regime are now open for innovation. 


Introduction 

Metamaterial research in the past decade offered an 
entirely new route to further enhance our capability to 
engineer material properties at will. Here, metamateri¬ 
als are artificially fabricated structures (often periodic, 
i.e., crystalline) which are designed so that they ex¬ 
hibit wave properties not observed with common mate¬ 
rials, e.g ., they can, in theory, bend electromagnetic [I], 
acoustic [2J, and even surface gravity waves [3] so as to 
achieve sub-wavelength focusing [4], create cloaks 
and attain shielding [7] . Other unexpected properties in¬ 
clude artificial magnetism [5j, negative permeability [2], 
negative refraction index [10], and hyperbolic dispersion 
m , to name a few. Such materials have allowed us to 
gain unprecedented control over a range of electromag¬ 
netic/optical and acoustic wave phenomena. In many 
ways metamaterials parallel the development of photonic 
and phononic crystals (optical and acoustic analogues of 
semiconductors) which also rely on small-scale structures 
for their properties. However, the major difference lies 
in the sub-wavelength nature of metamaterial structure. 
This enables us to summarize their properties in terms 
of permittivity and permeability (e, p) for electromag¬ 
netic waves, or bulk modulus and mass density (k, p) for 
acoustic waves, just as we would for any naturally oc¬ 
curring material. This is an enormous simplification for 
the design process, and research is now focusing on the 
realization of a new generation of metadevices [12j with 
novel and useful functionalities achieved by the struc¬ 
turing of functional matter on the sub-wavelength scale. 
Novel devices such as superlens |T5], hyperlens m , in¬ 
visibility cloaks [T51 fTB] . and plasmonic waveguides ca 
have been fabricated and tested in the past few years. 
The technology behind such metadevices is fairly well es¬ 
tablished in the low-frequency regime where inclusions 


have sizes much smaller than the wavelength of opera¬ 
tion. At these relatively low frequencies this is commonly 
obtained by assuming only monopole and/or dipole in¬ 
teractions, e.g., by utilizing conducting materials shaped 
as dipoles [T5J and split-ring resonators [19] , The ex¬ 
istence of resonances poses a considerable challenge to 
classical effective medium theories. This is because their 
basic principle is to minimize the scattering in the quasi¬ 
static limit, while the local resonances occur most often 
at longer wavelengths. 

Here, following these concepts, we develop and anal¬ 
yse a supra ^classical dynamic model of metamaterial 
response. There is an abundance of miscellaneous ef¬ 
fective medium theories [IH3J I23H26] , some quite recent 
mmm-, many of these works claim to be valid not 
only in the quasi-static limit but also at finite frequencies 
beyond the long-wavelength limit: a situation that hap¬ 
pens when the wavelength A is long in the host medium, 
while the wavelength in the particles, A 0 , can be small. 
(This is in contrast to the quasi-static limit where both 
A and Ao should be much larger than the size of the par¬ 
ticles.) Such extension to finite frequencies is sometimes 
denoted as the dynamic effective medium theory. How¬ 
ever, even this dynamic approximation relies exclusively 
on the monopolar and dipolar response of the scattering 
objects, which implicitly assumes long wavelengths. In 
this paper, this restriction is relaxed and the full effect 
of the ensemble of particles that constitute the effective 
medium is included, as higher diffraction orders are en¬ 
compassed. This will allow the design of new rnetade- 
vices working over a wider wavelength range. We shall 
illustrate this by solving a simple scalar problem in two 
dimensions, having applications not only in electromag¬ 
netics but also in acoustics and elasticity. The similarities 
between the equations of acoustics, elasticity and electro- 


* Dedicated with gratitude to the memory of Prof Yves C Angel 


1 This word comes from Latin and means above or beyond the 
limits of. 
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magnetics allow us to use some of the same techniques 
to solve problems in these seemingly disparate fields. 


View on classical homogenization 

The theoretical approach to the field of metamaterials 
is provided through dynamic homogenization techniques 
which relate the microstructure of a composite to the fre¬ 
quency dependence of its effective properties. The ma¬ 
jority of research interest in the area of metamaterials is 
restricted to periodic microstructures HU EH ( as the ar¬ 
rangement of molecules according to solid-state physics) 
which admit Bloch (or Floquet) waves as solutions and 
many different numerical algorithms have been developed 
(see, e.g., E2JI3Ej) for calculating the dispersive proper¬ 
ties of these waves. A popular route to determining these 
parameters is by the use of retrieval methods [34] [35] 
where the assumption is that local effective properties 
may be used to define periodic composites. The retrieval 
method leads to the refractive index n and the wave 
impedance Z, which defines the reflectivity of a semi¬ 
infinite slab. However, while simple in principle, such 
retrieval methods are limited to ordered arrays and often 
produce ambiguous results due to oversimplified initial 
assumptions of the bulk model [36) . 

Certainly engineers like structures and designs that fol¬ 
low some type of order. However, materials may be also 
amorphous and isotropic, and natural materials on the 
macroscopic level are quite often random in essence. It 
may well be that a random placement of complex parti¬ 
cles would be enough to produce emergent properties in 
the overall wave response and therefore give us a sample 
of metamaterial m ■ The effective behaviour of metama¬ 
terials whose microstructure is random depends strongly 
on the governing statistics of the random distribution. 
Effective properties may be determined by using the self- 
consistent effective medium methods for which a substan¬ 
tial body of literature may be found. Although variants 
exist, these methods often consider the scattering prob¬ 
lem of a coated particle embedded in a matrix which 
has the properties of the effective media. These prop¬ 
erties are then determined by requiring the vanishing of 
the effective forward-scattering amplitude /^ = /g£g 
and as such are formally restricted to the low-frequency 
and long-wavelength ranges. For examples where this 
method has been applied to electromagnetic, acoustic 
and elastic waves, see PH. Although the above self- 
consistent condition (/q^ = 0) is physically sufficient to 
describe the effective medium, two effective properties, 
i.e. ( e e ff,Heff ), cannot be determined ‘simultaneously 
and uniquely’ from the single condition. A supplemen¬ 
tary condition is needed; this prevents the application 
of effective medium methods to finding dynamic effective 
properties. Note however that the above condition is suf¬ 
ficient for wave propagation in a metamaterial in which a 
single material constant is involved, e.g. in dielectric me¬ 
dia. Another deficiency of many current enhancements 


of the effective medium methods is their failure to de¬ 
scribe the influence of the spatial distribution of parti¬ 
cles on the effective constitutive parameters. Such a de¬ 
scription is possible in the framework of a self-consistent 
scheme called the effective field method m and our work 
is within the framework of this scheme. One of the prin¬ 
cipal results of the effective field approach was an ade¬ 
quate definition of the coherent wave and a proof that it 
obeys a wave equation, i.e., a proof that, under certain 
conditions, a random distribution of scatterers can, for 
this purpose, be represented by an effective medium [10] . 
Most calculations proceed by assuming the existence of 
such an effective medium equation. 

The subject of the present work is the macroscopic 
dynamic behaviour of the above composite medium, i.e., 
random distribution of particles. More precisely, we shall 
describe a heuristic scheme for evaluating the effective 
properties of nretamaterials. The approach is based on 
the idea that a certain effective field acts on each particle, 
as a consequence of the presence of the other particles; 
hence, the name effective field method. The framework 
we develop is based on the Fikioris-Waterman sh sa 
and Waterman-Pedersen [8] formalism to evaluate the co¬ 
herent wave motion on both sides of a semi-infinite array 
of particles. 2 More specifically, we consider an averaged 
wave motion, where all possible configurations of parti¬ 
cles are weighted by appropriate pair-correlation func¬ 
tions. In contrast to the effective medium methods, we 
derive a fully dynamic model for the effective constitutive 
parameters, which retains all the relevant information 
(particle geometry and physical parameters) provided by 
an expanded multipole solution. As a result, the the¬ 
ory discussed in the following is more complete and po¬ 
tentially more useful than previous approaches to derive 
effective material parameters. 


Results 

Here, we consider two specific polarizations in electro¬ 
magnetism, transverse electric (te) and transverse mag¬ 
netic (tm). In addition, parallel to the electromagnetic 
example is the mathematically identical case of acoustics 
and anti-plane elasticity. 

We then consider these as two dimensional problems. 
Indeed, exploiting the physics common to many types 


2 The later reference, is the earliest work to our knowledge 
to predict explicit relations for the effective bulk parameters 
\£eff,neff) in the dynamic range. The authors have also pre¬ 
dicted negative frequency-dependent p. e ff a ' single-particle res¬ 
onances although the plots only displayed the positive values. In 
fact, they only noted that “the effective parameters vanish or 
diverge at certain frequencies” without further comment, which 
suggests that the results were considered curious at that time. 
Currently it is common to have negative effective parameters, 
and much research on metamaterials is focused on this area. 
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TABLE I: Relationships among electromagnetic, acoustic and 
elastic material parameters. 



Electromagnet ics “ 

Acoustics 

Elasticity 

m 

£ tm ^te 

P F 

1/G 

d 

^TM £ TE 

1/k 

P SH 


“Observe that the permittivity and permeability for a specific po¬ 
larization can be related to a pair of acoustic and elastic constants. 
For instance, (e TE ,/r TE ) «-> (p SH ,p p ) <->• (1 /k, 1/G). 


of wave propagation, the idea of metamaterials has been 
implemented successfully for acoustic and elastic waves. 
Many of the conclusions drawn from photonics research 
directly apply to acoustic waves and acoustic metama¬ 
terials due to the essential similarity of the governing 
equations in the two cases. Realizing analogous results 
for elastic metamaterials is complicated by the fact that 
the governing equation for elasticity admits both longi¬ 
tudinal and shear wave solutions which are capable of 
exchanging energy between each other. However, anti¬ 
plane elasticity is a special state of deformation with just 
a single non-zero displacement Held, similar to transverse 
electromagnetics. The governing equation common to 
electromagnetics, acoustics, and linear anti-plane elastic¬ 
ity is detailed in Appendix A. 

Effective constitutive parameters depend on many fac¬ 
tors including the intrinsic properties of the particles and 
the host matrix, their shape and topology. The latter de¬ 
termines how the particles are distributed in the matrix. 
The system considered in our study is composed of two 
isotropic phases: cylindrical particles of arbitrary shape 
randomly distributed in a host medium with propaga¬ 
tion constant k = md for some (possibly complex) 
parameters (m,d) of the medium. Depending on the ap¬ 
plication, these material parameters could be, e.g ., com¬ 
pliance (1/G) and density (p) for shear horizontal po¬ 
larized elastic waves or permittivity (e) and permeability 
(p) in electromagnetism; a number of useful relationships 
among these parameters are summarized in Table I. 

In Appendix B, we briefly review the effective field 
method. Subject to the quasi-crystalline approxima¬ 
tion, two equations are obtained for which the effective 
wavenumber K, of some coherent wave motion (either 
electromagnetic, acoustic, or elastic) and the effective 
impedance Z, are given in implicit form. Note that, 
whereas the dispersion relation for K, is polarization- 
independent, the effective impedance depends on the 
type of the incident wave. These equations are the start¬ 
ing point of all further developments. Observe that the 
particles have a size distribution and their relative posi¬ 
tions are described by an arbitrary cross-pair distribution 
function gij. Also, the size distribution is represented by 
Vj = v( a j)'i here a,j is the radius of the circular surface 
circumscribing a particle, and r) is the number of particles 
per unit area. 

Without loss of generality, we next assume the parti¬ 


cles are identical and have equal sizes aj = a. Here, we 
refer only to the final explicit solutions for the effective 
parameters m e ff and d e ff, which are expressed elegantly 
as 




d f- i+3 'i +3 4 + °(« 3 ) 


( 2 ) 


where e = 4 ^ 77 . The scalar coefficients (m 1 ,m 2 ) and 
(di, d 2 ) are given in matrix notation by 


and 


mi = e 4 Qe e 4 JQe, 

m 2 = e 4 QRQe e 4 JQRQe 

\2 / \ 2 


1 

4 k 2 V 


(e 4 Qe) — (e 4 JQe)' 


(3a) 


(3b) 


di = e 4 Qe + e 4 JQe, 

d 2 = e 4 QRQe + e 4 JQRQe 

,2 / \2 


1 

4 k 2 Y 


(e‘Qe) - (e 4 JQe)^ 


(4a) 


(4b) 


One can easily check that these equations are com¬ 
patible when m e ffd e fftC 2 = w 2 . Incidentally, we obtain 
K. 2 ~ k 2 + ee 4 Qe + e 2 e 4 QRQe + O (e 3 ), which is, as ex¬ 
pected, the second order expansion in e of the implicit 
wavenumber equation ( |B1 ). Note that all notations ap¬ 
pearing in Eqs. © and Q are introduced in the Ap¬ 
pendix. 

Results in classical multiple scattering theories are usu¬ 
ally defined in terms of the angular shape function fg for 
scattering of a plane wave by a single particle. It is useful 
to render yet another form of the coefficients ([3]) and Q 
in terms of fg. This is done by considering the line-like 
approximation: in addition to ea 2 < 1, we also require 
1. To render the results more tractable, the spa¬ 
tial distribution of particles is assumed to be isotropic 
and homogeneous, for which g,j (r) = g (r) = H(r — b). 
This describes a non-overlapping condition; here, g de¬ 
notes a pair-correlation function, H is the Heaviside unit 
function, and b = 2a is the diameter of the particles. Re¬ 
taining only the leading order term in kb of the multiple 


scattering matrix R, and using the definition (C2) for fg , 
we obtain 

e f QRQe = — —^R 0 and e 4 JQRQe = — (5) 


4 k 2 ' 


with 


JC„ = -/ d0cot(0/2)^ and 3 e a = fgf a -g. (6) 
77 Jo 
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By means of these approximations, we can infer the 
following closed-form constitutive relations 


m eff 

m 


1 + 2 ^ 
^[(sj-s8)-(^o-5c ff )], 


(7) 


if + + 

+ ^j[(s;-sS)-(«o+«.)]. (8) 

Apart from their dependence on k and e (or 77 ), the 
effective dynamic parameters (m e //,d e //) given by Eqs. 
([7]) and Q are all completely determined when the an¬ 
gular shape function fg for an isolated particle is known. 
If this scattering amplitude can be determined either an¬ 
alytically, numerically, or experimentally, then the effec¬ 
tive medium equivalent to the artificial composite is fully 
described. 

It is noteworthy that if one wants to study the be¬ 
haviour of effective parameters at high concentrations 
(where such expansions may not be valid) the general im¬ 
plicit equations detailed in Appendix B should be used 
and/or more accurate pair correlation functions should 
be considered. Neither incident wave nor boundary con¬ 
ditions have entered yet in the above description. Conse¬ 
quently, the results admit several solutions corresponding 
to different polarization states. In Ref. El (section SI), 
the expansions ([T]) and <[2j) [or ([7| and Q] are special¬ 
ized to electromagnetic, acoustic, and elastic scattering 
for long wavelengths (a -C A). This provides an addi¬ 
tional check on the correctness of the results obtained in 
this paper. A further check on the consistency of our 
method is provided in Ref. [44] (section S2). It is shown 
that the quasi-crystalline approximation is self-consistent 
and identical to coherent potential approximation [45] at 
least to second order in concentration, provided the ef¬ 
fective parameters are identified as those derived in this 
section. 


Discussion 


While the limiting cases considered in Ref. [33] (sec¬ 
tion SI) perform a check of the theory we have presented, 
they neglect some important features of the effective field 
method. Therefore, we address this problem numerically 
in order to illustrate the dynamic behaviour of the effec¬ 
tive parameters. In the following, the effective parame¬ 
ters (m e ff, d e ff) are calculated by using Eqs. @ and ( 2 ]), 
together with the Percus-Yevick pair-correlation function 
for hard disks [351 . 



FIG. 1: Sketch of the scattering of an electromagnetic wave 
by a circular cluster of cylindrical particles. 


Example illustration 

We consider a fibre bundle (or circular cluster of dielec¬ 
tric fibres with /to = 1 ) of effective radius r e // in vacuum. 
A plane electromagnetic wave is incident on the fibre bun¬ 
dle. A sketch is shown in Figure 1. There are 68 circular 
fibres each of radius a randomly distributed in the clus¬ 
ter and their volume fraction is 10.88%. The refractive 
index of the fibres is no(= \/^o) = 1-33 + O.Oli. Exact 
multiple scattering simulations 3 are compared with the 
effective medium model (i.e. equivalent homogeneous 
magneto-dielectric inclusion 4 with effective parameters 
(£ e ff,/j, e ff)). The multiple scattering results are aver¬ 
aged over different realizations of the fibres locations. 
With 500 realizations, the maximum error between the 
numerical model and the effective medium results is less 
that 0.5%, for the two cases illustrated. 

Figures 2 and 3 show the spatial maps of the near- 
field electric field amplitude \E Z \ for two different incident 
wavelengths, A = 2r e ///5 = 10a and A = 4r e ///25 = 4a, 
respectively. Figure 2 illustrates the response of the co¬ 
herent wave regarding the topology of the fibre-bundle. 
As expected, the waves are insensitive to the relative lo¬ 
cations of the fibres for long wavelengths. This is not 
the same for shorter wavelengths. A comparison of the 
results in Figure 3 indicates the agreement is excellent 


3 The analytical solution to Maxwell equations for scattering by 
multiple parallel cylinders has been described, e.g., in Ref. 09]. 

4 Note that, as expected, our results also predict an effective mag¬ 
netic permeability fj, e f f at finite frequencies (different from that 
in vacuum) in a system in which both the matrix and the parti¬ 
cles are non-magnetic. 
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Anisotropic metamaterials 



FIG. 2: Spatial maps of near-field electric field \E Z \ as a te 
wave is incident from the left. Exact multiple scattering sim¬ 
ulations for a single realization of fibres locations: left panels 
(A, C) - random array; right panels (B, D) - regular ar¬ 
ray. Top panels (A, B) - effective radius r e ff of the cluster 
is such that A = 2r e ///5 = 10a; bottom panels (C, D) - 
A = 4r e// /25 = 4a. 



o 



FIG. 3: Spatial maps of near-field electric field \E Z \ as a te 
wave is incident from the left. Left panels (A, C) - average 
over 500 different realizations (of exact multiple scattering 
simulations); right panels (B, D) - equivalent homogeneous 
inclusion (single scattering result) with dynamic effective pa¬ 
rameters. Top panels (A, B) - effective radius r e ff of the 
cluster is such that A = 2r e ///5 = 10a; bottom panels (C, D) 
- A = 4r e// /25 = 4a. 


even for the high frequency case considered (A = 4a). It 
is particularly encouraging that the agreement is excel¬ 
lent even inside the circular cluster. Observe that the a 
regular arrangement of fibres produces a result that is 
closer to the effective cluster for long wavelengths, than 
is the result obtained with a random realization of the 
fibres locations. We should note that although the com¬ 
parison in Fig. 3 is excellent, it may not always be so for 
other geometries of the fibre bundle. In a final section we 
detail various limitations and assumptions of our model 
and discuss other similar problems obtained previously. 


It is of considerable interest to discuss the possibility of 
realizing anisotropic metamaterials, that is, the material 
parameters are not scalars but tensors, with their prin¬ 
ciple components taking different values. Different from 
the anisotropy property of the material itself, we shall 
examine anisotropy originating from geometric asymme¬ 
try and consider a random array of elliptic cylinders of 
material parameters (mo,do). The x- and y -axes are set 
in the directions of the semi-minor and semi-major axes 
of the elliptic cylinders, with respective radii a x and a y . 
Due to the geometric arrangement of the elliptic cylin¬ 
ders and the symmetry of the scattered fields, the x- and 
y-directions can therefore be seen as effective principal 
directions. Proceeding essentially as detailed in Ref. [H] 
(section SI), we obtain, in the quasi-static limit 

+ ( 9) 

~ 1 + 2<$>M y +2<\> 2 (10) 

+ <t>V, (11) 

where </> = 7 Tr]a x a y is the volume fraction of the elliptical 

cylinders. The coefficients D and (M x ,M y ) are given by 


T> = — 1 and 

d 


M, 

M, 


1 (m 0 - m) (a x + a y ) 

2 m ofc)+m(M 


( 12 ) 


Observe that if a x = a y (= a, i.e., circular cross sec¬ 
tion) then M x =M. y (=A4, see Ref. [H], section SI), 
and m e ff, x = m eff,y The results of Eqs. ([£])-([TT]) show 
that only the effective property m e ff is a tensor with 
principal components m e ff, x an d m e ff,y, whereas d e ff 
is a scalar. This is consistent with results obtained re¬ 
cently in Ref. [6], for electromagnetic waves in the quasi¬ 
static limit. However, these results should be consumed 
with prudence. We show that, in general, both If t e ff 

and d e ff are tensors for arbitrary frequency and wave¬ 
length. To see this more clearly, let us consider the 
scattering of a TM wave by perfect electric conductive 
elliptic cylinders in vacuum 5 . From Table I, we infer 
that (m,d) corresponds to for tm waves; appro¬ 

priate identifications the resulting effective medium are 
implied. Figure 4 shows the effective permittivity ‘'e* e ff 
and permeability y, e ff tensors. Only the real part of 
these parameters is presented for brevity. The volume 
fraction cj) is fixed and equal to 67r%. It should be noted 
that the actual concentration </> = nr]a x a y cannot ex¬ 
ceed a x /a y in order to be consistent with our model, so 


5 This is equivalent to solving the Neumann boundary condition. 
In anti-plane elasticity, this condition corresponds to a cylindrical 
cavity with stress-free surface. 
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FIG. 4: Effective permittivity ^ e ff/e (= + ie^//) (top 

panels, A-C) and permeability *Jjteff /M (= Pe//+iMe//) (bot¬ 
tom panels, D-F) versus the wavelength for a TM-polarized 
field incident upon a random distribution of elliptical per¬ 
fect electrical conductive cylinders, for various aspect ratios 
a x /a y . Left panels (A, D) - a x /a y = 1; middle panels (B, E) 
- a x /a y = 0.75; and right panels (C, F) - a x /a v = 0.5. 

that (j> = ■nrja 1 < 1 when a x = a y (= a). The figure 
is intended to illustrate the variations of ^ e ff/e and 
*7i e ///£t as the wavelength A fa varies on the horizontal 
axis, for several aspect ratios a x /a y ; here, a = a x a y 
is the geometric mean of the semi-minor and semi-major 
axes, a x and a y . Observe that in the quasi-static limit, 
for A > 10a, where currently available model will be 
adequate, the principal components of yu. e // are visi¬ 
bly equal, i.e., p, e ff,x — l^eff^i regardless of the ratio 
a x /a y . This is as expected, given Eq. ( |TT| ). It is in¬ 
teresting that for shorter wavelengths (A < 10a), fi e ff,x 
and Heff,y become increasingly distinct as the ratio a x /a y 
decreases from 1 to 0.5, an effect not predicted by the ex¬ 
isting literature. This suggests a new route to the design 
of metamaterials with controllable anisotropic effective 
properties. 


Conditions of applicability 

The results in Figures 2 and 3 support the reliabil¬ 
ity of the effective material parameters resulting form 
the supra-classical dynamic homogenization procedure 
reported here. Note however that, although not appar¬ 
ent in the results, there is an approximation involved 
in replacing a finite-size heterogeneous composite with 
its homogenized equivalent, in addition to the reliability 
of the homogenization procedure (which ignores transi¬ 
tion region complications at the interface \x\ < a HD)- 
In practical terms, it means that for a finite sample of 
the random composite the applicability of dynamic ho¬ 
mogenization not only depends upon the frequency under 
consideration but also upon the phase of the composite 
at the boundary of the sample. An effort to quantify 


such an approximation is described in Ref. m ■ Here, 
the approximation results from truncating interfaces of 
a finite (or semi-infinite) 1-D periodic composite, the 
later being replaced with what are essentially its effective 
dynamic properties in the infinite Bloch-wave domain. 
Other questions will need to be answered relating to the 
shape and size of the scattering boundary, the effect of 
increasing the number of particles, and how many real¬ 
izations are required to determine both the near- and far- 
fields accurately. It is expected that as the bounded area 
increases, so does the uncertainty of the calculated field. 
An investigation in this direction is beyond the scope of 
the paper, however we refer the reader to the comprehen¬ 
sive numerical analysis (based on the quasi-crystalline 
approximation) reported in Refs. [HI [53] ■ The later 
references should come with a warning, as their analysis 
contains the implicit (and incorrect [54]) assumption that 
the dielectric permittivity is the only quantity of interest. 

Finally, let us note that the effective material pa¬ 
rameters derived in here (which are tensor values for 
anisotropic media) are not necessarily tied to the physical 
material parameters of any of the individual elements of 
the metamaterial. A rather critical survey discussing the 
link and the difference between these two concepts {i.e. 
effective versus characteristic material parameters), par¬ 
ticularly for the case of Maxwell’s equations, is presented 
in Ref. [55]. As evidenced in this survey and references 
therein, homogenization theories continue to attract at¬ 
tention and even controversy. It appears, from consid¬ 
ering exact reflection coefficients at oblique incidence (if 
one assumes that Fresnel-like formulae are always valid), 
that any effective material parameters that can be intro¬ 
duced in any theory would depend on the angle of inci¬ 
dence; broadly speaking, they would depend on the type 
of illumination. This means that these effective proper¬ 
ties do not necessarily relate solely to the bulk properties 
of the material itself; they can involve the material and 
the type of illumination. Relevant considerations in this 
direction are presented for periodic composites in Refs. 
[55][57]. A retrieval method extended to the arbitrary ori¬ 
entations of the principal axes of anisotropy and oblique 
incidence was presented in Ref. [58]. A discussion regard¬ 
ing modelling of the coherent wave propagation from the 
knowledge of the material properties along the principal 
axes only is elaborated in Ref. m- 

To summarize, a self-consistent multiple scattering ap¬ 
proach, which enables the dynamic homogenization of 
metamaterials in two dimensions is developed. The 
quasi-crystalline approximation is employed to break the 
hierarchy of increasing conditional probability densities, 
but otherwise the treatment is exact. In particular, the 
effective wavenumber and the effective impedance is ob¬ 
tained. These characteristics can then be used to de¬ 
termine the effective constitutive parameters of the ho¬ 
mogenised material. Whether the resulting effective pa¬ 
rameters represent a true bulk property of the metama¬ 
terial in the dynamic range is yet to be determined. The 
two natural approximations - dilute media and low fre- 
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quency approximations - show consistency, and, more¬ 
over, the quasi-static limit gives results reminiscent of 
the laws of Maxwell Garnett [2], Ament m , and Kuster 
and Toksoz m, respectively for electromagnetic, acous¬ 
tic, and elastic material parameters (see Ref. [33], sec¬ 
tion SI, for more details). The entire analysis described 
in this work is germane for alternative analytical proce¬ 
dures based on other scattering operators Q for an iso¬ 
lated particle. As shown in Ref. ( 31 ] (section S2), a fully 
self-consistent procedure may be based on a new kind of 
isolated scatterer problem. We have shown that the co¬ 
herent potential approximation, used in many previous 
works, is only an approximation of this procedure to the 
first order in the concentration of particles. 

The theory provided here offers exciting opportunities 
for researchers in different communities, ranging from 
seismic waves to the entire field of ultrasound research, 
and spanning radio frequency and optical engineering. 
In particular, metamaterial modelling in optics, physi¬ 
cal acoustics, and condensed matter physics may benefit 
from a rigorous, compact model for estimating more ac¬ 
curate and anisotropic effective medium parameters that 
homogenize artificial media. 
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Appendix 

A. Governing equation 


There are many examples of wave equations in the 
physical sciences, characterized by oscillating solutions 
that propagate through space and time while, in loss¬ 
less media, conserving energy. Examples include the 
scalar wave equation (e.g., pressure waves in a gas), 
Maxwell’s equations (electromagnetism), Schrodinger’s 
equation (quantum mechanics), elastic vibrations, and 
so on. From a mathematical viewpoint, all of these 
share certain common features. In the following, we shall 
briefly identify the similarities between three types of 
such waves, in two dimensions: electromagnetic waves, 
anti-plane elastic waves, and acoustic waves. Electro¬ 
magnetic waves are quite different from acoustic and elas¬ 
tic waves in that they can travel through vacuum. How¬ 
ever, from an algebraic perspective, all three types of 
waves can be described by a unique scalar equation and 
hence these disparate phenomena can be studied simul¬ 
taneously. The prototype problem consists of the wave 
equation at fixed angular frequency to 


V- (m 1 (r) (r)) + w 2 d (r) %p (r) = 0, (Al) 


where c = l/V md is the phase velocity of the wave for 
some parameters m (r) and d (r) of the medium. In the 
above, the scalar wavefunction ip (r) e _1 “* corresponds to 
some physical field. 

The problem considered here is reduced to points in the 
x-y plane (i.e. the cross section plane of our scattering 
geometry), which in polar coordinates are r = (r, 9); here, 
6 is measured from the positive x-axis. Let us first con¬ 
sider the two important modes for electromagnetic wave 
propagation: the transverse electric field and the trans¬ 
verse magnetic field. These modes are closely analogous 
to anti-plane shear in elastodynamics and to acoustic 
waves. Let us assume the medium is isotropic and has di¬ 
electric permittivity e and magnetic permeability p that 
are independent of z. A transverse magnetic (tm) field 
is a special solution of the Maxwell’s equations that has 
the form H = ip (x, y) L, and the electrical field iuie'E = 
— V x H = — ^i y . A transverse electric (te) field is 

another special solution of the Maxwell’s equations that 
has the form E = ip (x, y) i z , whereby the magnetic field 
is given by iw/iH = V x E = y — ^i x . Observe that 
simple knowledge of the scalar wavefunction ip suffices for 
the determination of the x and y components of the elec¬ 
tric and magnetic fields, for the two polarizations. Hence, 
Eq. (Al) is the governing wave equation for electromag¬ 


netic waves provided that (m,d,V0 ^ ( e ; A Hz) for TM 
waves, and (m,d,z/>) «->• (p,s,E z ) for TE waves. The gen¬ 
eral solution independent of z is a superposition of the 
TE and TM solutions. This can be seen by observing that 
the Maxwell’s equations decouple under this condition 
and a general solution can be written as ( H x , H y , H z ) = 
(H x , H y , 0) + (0, 0, H z ), where the second term represents 
the tm solution. The first term is of course the te so¬ 
lution because Vx (H x , H y , 0) = (0,0 ,^H X - -§^H y ), 
which implies E x = E y = 0 as expected. Let us now 
consider the case of anti-plane shear strain which is a 
special state of deformation where the displacement field 
is given by u = ip (x, y) i z . This is an out-of-plane mode of 
deformation and is analogous to transverse electromag¬ 
netic wave propagation. In the linear regime, an isotropic 
elastic medium is characterized by its density p and the 
Lame elastic constants G = p and A: G is the shear 
modulus (notation used to distinguish from the perme¬ 
ability p employed in electromagnetics) and A = n ~ G 
where k is the two dimensional bulk modulus. Hence, Eq. 


(Al) is the wave equation for anti-plane shear provided 


that (m,d,^>) (1 /G,p,u z ). Shear waves that satisfy 
this equation are also called SH (shear horizontal) waves, 
particularly in seismology. Finally, let us consider the 
acoustic wave propagation in an isotropic medium. For 
an inviscid fluid or gas, the shear modulus G is zero - 
and A is just the bulk modulus. In this case, replacing 
ip with the pressure p = — • u, we obtain precisely 

the acoustic wave equation (Al I for p (or compressional) 
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waves if (m, d) (p, 1 /k). 

In essence, the solutions of the three problems con¬ 
sidered will lead to similar conclusions if we make the 
appropriate interpretation of the quantities involved. 


B. Effective field method 

Suppose that discrete particles of cylindrical geometry 
are randomly and uniformly distributed in a half-space 
defined by {x > 0}. The particles need not be circular, 
provided that each of them can be contained in a circum¬ 
scribing circular surface of radius aj (with an axis of rev¬ 
olution parallel to the z-axis); their number density is ijj. 
Both the particles and the matrix are made of isotropic 
materials. Let a plane wave ip = exp [i (kx — uot)] of unit 
amplitude propagate with wavenumber k in the matrix 
along the ^-direction. When this wave propagates in the 
composite material, multiple scattering occurs. Either 
propagation or diffusion, or a combination of the two 
phenomena is observed, depending on the frequency as 
well as on the geometrical and material properties of the 
composite. Assuming that propagation occurs, one can 
describe the coherent wave motion in the composite by a 
complex-valued wavenumber 1C. The fundamental equa¬ 
tion for configurational averages of the exciting and total 
fields for scalar wavefunctions has been derived in detail 
in Refs. H3 ED ES ED]. The quasi-crystalline approx¬ 
imation |6lJ is used to truncate the hierarchy of equa¬ 
tions (Foldy-Lax hierarchy) so that only the correlation 
between every two particles is considered. We obtain the 
implicit dispersion equation for the effective wavenumber 
K. of the coherent wave exp [i (ICx — tot)}, 

1C 2 = k 2 + Yhjtj-F ( a j) , (Bl) 

where = 4 ^ 77 ^, and the effective scattering amplitude 
T is given by 

F{a 0 ) = e*(Q” 1 - E^R^r'e, (B2) 

e = (1,1,... Y is a constant unit vector. The shorthand 
notation = R (bij) and Qj = Q (aj) has been used. 
The infinite square matrices R ( bij) and Q (aj), have el¬ 
ements 


n nv (bij) = Vn ^ h _ l \l 1 l +M n - v (ICbij), (B3) 


where bij ^ + aj, and 


Qnv (aj) — . b nv T n (aj ). (B4) 

17T 

Here, <5 i; y denotes the Kronecker delta, and V# and Mi are 
given by 


'Pi 0) = y [ zH \ L> ( x) f~ z Ji (z) - xJi (z) £Hj 1J (x) 


rd) 


r(D 


17r 


Me (z) = ^ I dr [gij (r) - 1 }rH ( p (kr) J e (zr/hj) 


IJ ) 5 


where Ji and H are the cylindrical Bessel and Hankel 
functions, respectively, and x = kbij. The function is 
the cross-pair distribution function of two particle species 
(with sizes cii and aj), and satisfies the non-overlapping 
condition: g^ (r) = 0 for r < bij ; also, if the distance 
between particles tends to infinity, then the correlation 
between their locations disappears, i.e., lim gij (r) = 1. 

r—> 00 

The scattering coefficients T n (aj) in Eq. (B4) depend on 
frequency, size aj, as well as on the properties of the par¬ 
ticle and those of the matrix material; they are evaluated 
by imposing app ropriate boundary conditions at r = aj. 

Equation (Bl) follows directly from a Lorentz-Lorenz- 


type law, and is an exact expression for the effective 
wavenumber, subject to the quasi-crystalline approxima¬ 
tion. It is of interest to note how various physical aspects 
are embedded in this equation. The scattering matrix Q 
describes the response of a single particle to a plane in¬ 
cident harmonic wave with wavenumber k, and contains 
all the scattering behaviour in terms of particle geome¬ 
try and physical parameters. The effective wavenumber 
K, only appears in the matrix R, which is defined by the 
spatial arrangements of particles, and accounts for mul¬ 
tiple scattering. Should the distribution of particles be 
regular, the quasi-crystalline approximation is exact, in 
which case the multiple-scattering matrix R can be re¬ 
duced to a well known lattice sum. 

The theory described above is now complete insofar as 
behaviour within the medium is concerned. It is also of 
interest, however, to calculate the effective impedance Z, 
which defines the reflectivity of the half-space {x > 0} - a 
quantity which may be measured directly. Following the 
derivations in Refs. mm, the coherent reflected field 
(ip) = fHexp(— ikx) at the half-space boundary can be 
obtained explicitly, with the reflection coefficient defined 
as 


Sy e jF-n (aj) 

4 k 2 + tjMo (aj )' 


(B5) 


Here, represents the average (coherent) back-scattered 
amplitude at normal incidence in the domain {x < 0}. 
The effective scattering amplitudes, Jo and M, corre¬ 
spond to coherent waves scattered in the forward and 
backward directions, respectively, and are given by 


Jo (aj) = etQjVj and J' 7r (aj) = e*JQ jVj, (B6) 

where J = {5 nv cos n7r} is a diagonal infinite matrix. The 
infinite eigenvector v ; , associated with the wavenumber 
equation, follows from an Ewald-Oseen-type extinction 
theorem, with the result 








2k 

Vj ~ K + k 


(l ,Q, R,,) e, 


(B7) 


where I is a unit infinite matrix. Martin [62. has ob¬ 
tained a formula for for obliquely incident waves on a 
half-space of circular scatterers; it can be shown that at 
normal incidence the result in his Eq. (39) gives agree¬ 
ment with Eq. (B5). The behaviour of the fields across 
interfaces was also examined in Refs. 


It was 

found that the fields themselves are continuous but the 
slopes are discontinuous. Using the estimate for the slope 
discontinuity, effective constitutive parameters can be de¬ 
rived, as shown in Refs. [62J [Me Equation ( |B5| ) is an ex¬ 
act formula for the reflection coefficient 6 ; it can be used 
to determine effective parameters (m e //,d e //) uniquely. 
It is often assumed that the effective medium correspond¬ 
ing to the distribution of particles may be described as 
a homogeneous medium from the standpoint of coherent 
wave propagation - the homogenized equivalent having 
the effective dynamic properties of the composite. In the 
following, we shall use this analogy, whereby the reflec¬ 
tion coefficient at the interface between the homoge¬ 
neous medium and the homogenized equivalent, may be 
written (as is standard) in terms of impedances (result¬ 
ing in a Fresnel-like formula). 7 Then, equating the result 
with Eq. (B5), the effective impedance Z can be explic¬ 


itly calculated. As expected, the effective impedance is 
different for different polarizations. Two cases are possi¬ 
ble, with the following results: 


Z m = Z/z and Z d = z/Z, 


(B8) 


where (Z m ,Z d ) are defined in Eqs. (B8)-(B9). Observe 
that, by using the definition, 


/C = k+ (aj) 


(Bll) 


in Eq. (BIO), the resulting parameters (m e //,d e //) can 


be expressed explicitly in terms of the effective forward 
and back-scattering shape functions, Tq and J- n . 

To conclude this section, we consider the line-like ap¬ 
proximation of the constitutive parameters (BIO). For 
this, the size of the particles is assumed small 
to the incident wavelength (aj -C A). 


compared 
At leading or¬ 
der, the single-scattering operator Q is compact and has 
only three eigenvalues of finite size (related to terms with 
n = 0, ±1). Furthermore, the infinite multiple-scattering 
operator R is reduced to a rank 3 matrix. Omitting the 
details, we find for circular cylinders (with Ti = T_i), 


n Q7 

m 


fc2 + h Ej o T i (aj) 
fc2 - h o T i («t) ’ 




(B12) 

(B13) 


where 


It can be shown that in the quasi-static limit (aj <C Ao) 
the effective property m e ff is reminiscent of the laws 
of Maxwell Garnett m, Ament [12], and Kuster and 
Toksoz na, in two dimensions, respectively for elec¬ 
tromagnetic, acoustic, and elastic material parameters. 
(This is further described in Ref. [S3, section SI.) On 
the other hand, the effective property d e // reduces to the 
simple and inverse rules of mixtures, depending on the 
physical model under consideration, and as seen from Eq. 
(B13) is linear in ej. 


„ _ 4* 2 + Ej U Fo (a,-) - (%)] 

Z 4kZ + E j e j [J r o(a j )+^(a j )Y [ 

Here, z = -y/m/d is the impedance of the matrix; the 
superscripts ‘m’ and ‘d’ correspond to different physical 
situations, as we shall see below. We can now state our 
most general expressions for the effective dynamic con¬ 
stitutive parameters ( m e jf , d e //), 

= £Z m and ^ (BIO) 

m k ok 


6 Note however that during the derivation, complications in the 
transition region —a<x<a near and on both sides of the 
interface have been ignored m- 

7 For instance, for acoustic waves £R = (Z — z) / {Z + z)\ this re¬ 

sult implies the continuity of pressure and normal velocity at the 
interface; for anti-plane elastic waves = — (Z — z) / (Z + z)\ 
here, the continuity of the out-of-plane displacement and the 
corresponding stress are implicit. Similar results in electromag¬ 
netics are known as Fresnel relations (for TE and TM waves). 


C. Explicit second order approximations 

At low concentrations (cjaj <C 1), the dispersion equa¬ 
tion is explicit, and reduces to the well-known formula 

SOI 

/c 2 ~ fc 2 + Ej e jfo ( a j) ? (ci) 

where the forward-scattering amplitude fo is given by 
/o (aj) = e f Qj-e. More generally, the angular shape func¬ 
tion fg for each particle is defined, in terms of Fourier 
series, as 8 

f S (a J ) = -J2 T n( a j)e in0 - (C2) 

17r c ' 

n 


Note the shorthand notation = y] is used throughout. 
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An expansion of the dispersion equation (Bl) to the 
second order in concentration results in 


K," — fc 2 +£ i e J -e t Q j e+£ i jG Cjc f Q,R, v / Q,e +O (eiCjek), 

(C3) 

where the matrix = R(&y) = lim R )? and has ele- 

/C —y k 

ments 


where the back-scattering amplitude is given by 
f n ( a,j ) = e 4 JQje. The second order approximation is 
too long to warrant including here. For completeness, we 
also give the following results, in terms of Fourier series, 


e QiR'ijQje 



y Unu (hj) T n (a/) T v (ttj), 

n,v 


(C6) 


n. 


nv (hj) = Me, (x) + [ (L 2 - X 2 ) J e ( X) Hf> (x) 


-* 2 a h J ^(x)£Hj 1] (*)], 


(C4) 


with l = n — v, and x = kbij. Note that for spatially 
uncorrelated particles, AT? ( x ) = 0. 

For the effective impedance of Eq. (B91, at first order 
in concentration, we have 


2 = z - (Oj). 


(C5) 


«> .l Q R , Q e ^ (-l) n 1l nu ( bij) T n (cq) T v ( a 3 ) . 

n,v 

(C7) 

These expressions can be easily approximated in the 
low frequency limit by observing that, to leading or¬ 
der in (kbij), and for uncorrelated particles, TZ nl , = 
\n — v\ /2k 2 . The results obtained here have been used to 
derive the analytic formulae presented in the main text. 
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Supplementary text 


SI. Consistency check 

It is straightforward to expand the low-frequency for¬ 
mulae obtained in the Appendix. Instead, the expansions 
(1) and (2) [or (7) and (8)] are approximated for long 
wavelengths (a <C A). 

Small cylinders behave as a combination of a monopole 
(or a source) and a dipole: this is a generic situation, 
and the leading order contribution to the angular shape 
function fg only involves the scattering coefficients Tq (w) 
and Ti ( uj ), i.e. 


fe* A {T 0 + 2Ti)cos0. (SI) 

17r 

The integrals in Eq. (6) reduce to 

Wo = 4n (To + Ti) and (S2) 

7T Z 7r z 

Finally, using these results, the effective parameters 
may then be expressed as 


m eff 


leff 


~ 1 + e- —tttTi - e 2 


Si 1 + 67 


ink 2 

1 


2 fc 4 


rjiZ 
1 1 > 


i7rfc 2 


To 


(53) 

(54) 


which agree to 0 (e 2 ) with the general estimates (B12) 
and (B13) in Appendix B, as expected. We can see that 
d e ff and vn e ff are related to monopolar (n = 0) and 
dipolar (n = 1) scattering coefficients, respectively. Note 
that the results in Eqs. (fS3l) and (fS4l) need not corre¬ 


spond to the quasi-static limit, because wavelength (and 
1/e length, if losses are present) within the particles is 
as yet arbitrary relative to the particle size. Many previ¬ 
ous effective medium results |SlIiS6) correspond to such 
dynamic approximations. If we further assume that the 
wavelength within the particles is small (a <C Ao), it can 
be shown that to the leading order in (fca), 


where 


T 0 ~ jk 2 a 2 V and Ti ~ ^ k 2 a J M , 


(S5) 


d m 0 + m 


(S6) 


The quasi-static relations (S5) with coefficients (S6) 


hold for a variety of boundary value problems, equiv¬ 
alent across electromagnetics, acoustics, and elasticity. 
This includes the Neumann boundary condition by set¬ 
ting do = rno = 0. As is well known, the Dirichlet condi¬ 


tion is atypical, and special care is needed. 9 In electro¬ 
magnetism, Dirichlet or Neumann boundary conditions 
(depending on the polarization in question) describe in¬ 
clusions that are perfectly conducting. 

In terms of T> and M, Eqs. (S3) and (lS4| yield 


~ 1 + 2(j)M+2<t) 2 M 2 , (S7) 

m 

=* 1 + 4 >v, (S8) 

d 

where <j> = ea 2 /4 (= nrja 2 ) denotes the fractional vol¬ 
ume occupied by the particles. These equations reduce 
to different forms in the quasi-static limit depending on 
application. Let us now specialize them to electromag¬ 
netic, acoustic, and elastic scattering, in succession. 


Electromagnetic waves 

There is a vast literature describing the many ap¬ 
proaches to calculate effective-medium electromagnetic 
parameters [S8] . and many of the existing theories are 
closely related to models developed in the late 1800s and 
early 1900s. We note in particular the expressions for the 
effective permittivity £ e ff obtained in three dimensions 
by Maxwell Garnett [S9] and by Bruggeman |S10j that, in 
turn, are closely related to the older Lorentz-Lorenz for¬ 
mula for time-dependent electric fields and the Clausius- 
Mosotti equation for static fields m ■ The Bruggeman 
formula has the special property that it treats the par¬ 
ticles and the environment symmetrically. However, this 
results in a quadratic order in <fi that is different from 
the expansion of Maxwell Garnett rule. Conversely, our 
result for e e f f is consistent with the rule of Maxwell Gar¬ 
nett. Indeed, consider a mixture where small magneto¬ 
electric particles are embedded in a host environment of 
permittivity (permeability) e {fi ). The complex permit¬ 
tivity of the particles is £o and their permeability is /to- 
Then, for TM waves, we have 


e% f ~ e + 2 (j>e— —^+2</> 2 £ 
e)i £ 0 +£ 

Z 4 ™/ - M + </> Oo - fi) ■ 


£q-£ 
£o + £ 


2 


(S9) 

(S10) 


Should the incident wave be electric in nature, £ and /i 
would have been interchanged. In the language of the di¬ 
electric problem, Eq. (S91 is reminiscent of the Maxwell 
Garnett estimate (more precisely, a small-</> approxima¬ 
tion of the Maxwell Garnett rule in two dimensions). On 
the other hand, the effective permeability of Eq. (S10) is 
given by a simple rule of mixtures. 


9 A particular feature of the Dirichlet problem is the presence of 
log ka in the asymptotics of the solution. The interpretation of 
this effect depends on the physical model under consideration. 
For investigations in this direction, see, e.g., Ref. m- 
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Sound waves in a compressible fluid 

In the acoustics context, consider a fluid-particle mix¬ 
ture and let (p, n) and (po, kq) be mass densities and bulk 
moduli of their respective phases. The particle material 
can be, e.g. solid or fluid. Then, 


Plff 

K eff 


p + 2 <j>p— — ~+2(f) 2 p 
Po + P 


1-0 

K 



f Po- p V 
\Po+p) 


(511) 

(512) 


The mass density p v e ^ is analogue to that obtained 
in three dimensions by Ament m ( more precisely, a 
small-0 approximation of the Ament-estimate in two di¬ 
mensions). Equation (S12) is recognized as the Reuss 
average for the effective bulk modulus n e ff- 


below is not limited to only these geometries but appli¬ 
cable, in principle, to any other particle shapes. From 
Eq. (Bl), we consider the effective wavenumber at the 
second order in concentration as, 

A 2 = k 2 ff + ee*Q e //e + e 2 e l Q e ffR e ffQ e ffe + O (e 3 ) , 

(S15) 

where the effective operator Q e f f corresponds to a coated 
particle excited by the coherent motion exjp(ik e ffx). The 
multiple-scattering matrix R e // is given by Eq. (C4), 
but with k replaced by k e ff. The self-consistent scheme 
now assumes that /C = k e ff. From a physical point of 
view this means that the coherent wavefield the compos¬ 
ite medium coincides with the wavefield propagating in 
the effective medium. Hence, the medium can be con¬ 
sidered as homogenized since there is no scattering in 
the outer effective medium. This results in the following 
non-linear equation for k e f f , 


Anti-plane elastic waves in a solid composite 

Consider the elastodynamic problem of anti-plane 
shear scattering, and let the host and particles have shear 
moduli G and Go, respectively. Again, the particles can 
be made of, e.g. solid or fluid material. For solid parti¬ 
cles, we obtain 


G 

G^f 


n SH 

Peff 


1 i ‘TJ.G ~ Go , 0 ^2 

1 + 2 ^gTGo + ^ 
p + <t> (po - p) ■ 


(G-G 0 
\G + G 0 
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(513) 

(514) 


The effective shear modulus G e ff is analogue to that 
obtained in three dimensions by Kuster and Toksoz [S13j 
(more precisely, a small-0 approximation of the Kuster- 
Toksoz estimate in two dimensions). Observe that the 
effective mass density p™ j in the case of anti-plane elas¬ 
ticity is also given by a simple rule of mixtures. 


S2. Self-consistent effective field method 

Self-consistent methods for the problem of scalar wave 
propagation through a medium with many particles may 
be found in the works of Maxwell and Rayleigh. During 
more than a century, in a number of works, these methods 
were extensively developed and used for the solution of 
various wave propagation problems. The present general 
results may be recast in terms of a dual-layer scattering 
operator Q e f / corresponding to a coated particle embed¬ 
ded in an effective medium with the properties m e ff and 
d e ff. The coating is made of the original matrix material 
with the properties m and d. For simplicity, we assume 
the particles are circular cylinders and have equal sizes 
aj = a. The radius c of their coating defines the volume 
fraction to be 0 = a 2 /c 2 . Note that the theory presented 


e Qe// e — ee Qe//R-e//Qe//G + O (e 


(S16) 


It is thought that the use of the self-consistent scheme 
(S16) applied to the effective wavenumber (Bl) can im¬ 
prove the accuracy of the results while the concentration 
of particles increases. Note that since Q e // is a transcen¬ 
dental function of the unknown fc e //, explicit solutions 
can only be obtained subject to low-frequency approxi¬ 
mation. Writing the dual-layer forward-scattering ampli¬ 
tude JqG 0 f a coated particle, as JqG = e 4 Q e //e, Eq. 
(S16) is then reduced to 


f e ff ~ _ 
jo — : 


2e 


2 k 2 


^ r Tf f ( r T 1 ^f f 

— T i (V 

eff V 


-t: 


•eff 


(S17) 


The effective scattering coefficients TqG anc f can 
be calculated in terms of 0 (= ea 2 / 4), as 


rreff m , 2 2 

T ° -4 t kef f a 


] e// 


(1 + 02 ?) - 1 


T Gf~ — k : 


(S18) 


(S19) 


40 e G (1 + cj)M) m+ (1 - <f>M) m e ff ’ 


where M. and V are defined in Eq. (S6). We find directly 


that T-G f = G(0 3 ), if one uses the small-0 estimate 
(S7). More generally, without restriction on 0, we obtain 
TqG = = o, for the following effective parameters, 


m eff 

m 

d eff 

d 


1 + (pM 
1 -0A4’ 

1 + <f>D. 


(520) 

(521) 


Here, and d e // are the quasi-static limits of Eqs. 

(B12) and (B13). This result is significant, since (S17) 
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reduces to the coherent potential approximation of solid 
state physics, regularly used in electromagnetics, acous¬ 
tics and elastodynamics, and is equivalent to /q^ = 0 . 
The coherent potential approximation appears to be an 


approximation of the solution to Eq. (S17) to first order 


in concentration. This leads to two important conclu¬ 
sions. First, the application of the self-consistent scheme 
to the quasi-crystalline approximation does not change 
the result ( e.g. we obtain the same effective wavenumber 
in both cases, if and only if the effective parameters de¬ 
rived in this paper are employed), and second, it reduces 
exactly to the coherent potential approximation, at least 


to second order in concentration. The later is in contrast 
with the findings of Ref. ihh, where the wrong effective 
parameter was employed in the self-consistent scheme, 
that is rr\,,ff = m + </> (m 0 — m) (describing the effective 
mass density of a fluid-particle mixture). Finally, let us 
note another important aspect of the self-consistent effec¬ 
tive field method presented in this section. The new self- 


consistent scheme (S16) includes the influence of the spa¬ 


tial distribution of particles (through the matrix R e //), 
whereas the coherent potential approximation does not. 
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